#!/usr/bin/perl
use Cwd 'abs_path';
use Getopt::Std;
#Input Arguments
getopt("d:t:a:f:o", \%args);
#my $dir_or = $args{d};

my $tped= $args{t};
my $arraytype= $args{a};
my $tfam = $args{f};
my $outdir = $args{o};
my $dirtemp = $args{d};
chomp($dirtemp);
if($dirtemp eq "")
{
        die "entered temp dir name  is empty\n";

}
if (!(-e $tfam)) {
        die "input tfam file does not exist\n";
}
if (!(-e $tped)) {
        die "input tped file does not exist\n";
}
if($tped eq "")
{
        die "entered tped file name  is empty\n";

}
else
{
@tped_ar=split(/\//,$tped);
$tped = pop(@tped_ar);
$dir_or = join('/',@tped_ar);
}
if($tfam eq "")
{
	die "entered tfam file name  is empty\n";
}
if($outdir eq "")
{
	die "entered output directory name  is empty\n";
} 

#to know the path of this script to grab working directory location
$line = abs_path($0);
chomp $line;

@DR_array = split('/',$line);
pop(@DR_array);
$dir = join("/",@DR_array);
$round = sprintf("%.0f", rand()*time());
#$round =1218378123;
$rounded = "temp".$round;

#$rounded = "temp265747191";
#creating a temp directory if not presen
system("mkdir -p $dirtemp");
system("mkdir $dirtemp/$rounded");
system("cp $tfam  $dirtemp/$rounded");
@tfam_ar=split(/\//,$tfam);
$tfam = pop(@tfam_ar);
system("mv $dirtemp/$rounded/$tfam $dirtemp/$rounded/unprocessed_beagle_input.tfam");
system("cp $dir_or/$tped $dirtemp/$rounded/unprocessed_beagle_input.tped");
$ori_tped = $tped;
$tfam ="unprocessed_beagle_input.tfam";
$tped = "unprocessed_beagle_input.tped";
@array_arraytype = split(',',$arraytype);
open(BUFF,"$dirtemp/$rounded/totalsnps");
close(BUFF);
open(BUFF,"$dirtemp/$rounded/totalsamples");
close(BUFF);
#if(@array_arraytype > 1)
#{
	for($i=0;$i<@array_arraytype;$i++)
	{
		$refpath = "$dir/bin/FwdStrand_Naresh/ref/".$array_arraytype[$i]."/FWD_REF_INFO";
		print $refpath."\n";
		system("mkdir $dirtemp/$rounded/$array_arraytype[$i]");
		system("cut -f2 $refpath > $dirtemp/$rounded/$array_arraytype[$i]/snps_platform_$array_arraytype[$i]");
		system("$dir/bin/plink/plink --tfile $dirtemp/$rounded/unprocessed_beagle_input --extract $dirtemp/$rounded/$array_arraytype[$i]/snps_platform_$array_arraytype[$i] --make-bed --out $dirtemp/$rounded/$array_arraytype[$i]/temp");
		if($array_arraytype[$i] =~ m/AFFY/)
        	{
                	$geno = 0.1;
                	$mind = 0.1;
        	}
        	else
        	{
                	$geno = 0.05;
                	$mind = 0.05;
        	}
		print "geno: $geno mind: $mind\n";
		#running plink to get the missing information
        	system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/$array_arraytype[$i]/temp --missing --out $dirtemp/$rounded/$array_arraytype[$i]/miss");
        	if (!(-e "$dirtemp/$rounded/$array_arraytype[$i]/miss.imiss")) {
                	die "something wrong in the plink command miss \n";
        	}
		#Remove SNPs that fail in all the samples and Remove Samples with >15% missing. 
        	system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/$array_arraytype[$i]/temp --missing --geno 0.999 --mind 0.15 --make-bed --out $dirtemp/$rounded/$array_arraytype[$i]/intialqc");
        	if (!(-e "$dirtemp/$rounded/$array_arraytype[$i]/intialqc.fam")) {
                	die "something wrong in the plink command geno\n";
        	}
		
		#remove SNPs that failed in >=5% for illumina and >=10% for affy of remaining samples
        	system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/$array_arraytype[$i]/intialqc --missing --geno $geno --make-bed --out $dirtemp/$rounded/$array_arraytype[$i]/geno");
        	if (!(-e "$dirtemp/$rounded/$array_arraytype[$i]/intialqc.fam")) {
                	die "something wrong in the plink command geno\n";
        	}
		
		#Flag Samples with a lot of missing>(5% Affy and 2% Illumina).
        	system("$dir/bin/plink/plink --bfile $dirtemp/$rounded/$array_arraytype[$i]/geno --missing --mind $mind --make-bed --out $dirtemp/$rounded/$array_arraytype[$i]/mind");
        	if (!(-e "$dirtemp/$rounded/$array_arraytype[$i]/mind.fam")) {
                	die "something wrong in the plink command mind\n";
        	}
		
		#imputing the sex information to cross check with the given one	
		system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/$array_arraytype[$i]/mind --impute-sex  --recode --transpose --out $dirtemp/$rounded/$array_arraytype[$i]/imputesex");
        	if (!(-e "$dirtemp/$rounded/$array_arraytype[$i]/imputesex.tfam")) 
		{
                	die "something wrong in the plink command impute sex\n";
        	}
        	system("mv $dirtemp/$rounded/$array_arraytype[$i]/miss.imiss $dirtemp/$rounded/$array_arraytype[$i]/samplesmiss.qc");
        	#system("mv $dir/temp/$rounded/samples.qc $outdir/samples.qc");
        	#system("rm $dirtemp/$rounded/$array_arraytype[$i]/miss.*");
        	#system("rm $dirtemp/$rounded/$array_arraytype[$i]/geno.*");
        	#system("rm $dirtemp/$rounded/$array_arraytype[$i]/mind.*");
        	system("mv $dirtemp/$rounded/$array_arraytype[$i]/imputesex.sexcheck $dirtemp/$rounded/$array_arraytype[$i]/sexcheck.txt");
        	system("mv $dirtemp/$rounded/$array_arraytype[$i]/imputesex.tped $dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tped");
        	$sys = "perl $dir/bin/perl_replace_originalsex_predicted_sex.pl $dirtemp/$rounded/$array_arraytype[$i]/sexcheck.txt $dirtemp/$rounded/$array_arraytype[$i]/imputesex.tfam $dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tfam";
        	print $sys."\n";
        	system($sys);
        	#system("rm $dirtemp/$rounded/$array_arraytype[$i]/imputesex.*");
        	#system("cp $dirtemp/$rounded/unprocessed_beagle_input.tfam $dirtemp/$rounded/processed_beagle_input.tfam");
        	$sys = "perl $dir/bin/sexcheck.pl -i $dirtemp/$rounded/$array_arraytype[$i]/samplesmiss.qc -d $dirtemp/$rounded/$array_arraytype[$i]/sexcheck.txt -v $mind -o $dirtemp/$rounded/$array_arraytype[$i]/samples.qc";
        	print $sys."\n";
        	system($sys);
		$sys = "cut -f2 -d ' ' $dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tped >> $dirtemp/$rounded/totalsnps";
		print $sys."\n";
		system($sys);
		$sys = "cat $dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tfam >> $dirtemp/$rounded/totalsamples";
                print $sys."\n";
                system($sys);		
		$sys = "mv $dirtemp/$rounded/$array_arraytype[$i]/samples.qc $dirtemp/$rounded/$array_arraytype[$i].qc";
		$qc[$i] = "$array_arraytype[$i].qc";
                print $sys."\n";
                system($sys);

	}
	system("mv $dirtemp/$rounded/unprocessed_beagle_input.tped $dirtemp/$rounded/unprocessed_beagle_input_ori.tped");
	system("mv $dirtemp/$rounded/unprocessed_beagle_input.tfam $dirtemp/$rounded/unprocessed_beagle_input_ori.tfam");
	if(@array_arraytype == 1)
	{	$sys = "mv $dirtemp/$rounded/$array_arraytype[0]/unprocessed_beagle_input.tped $dirtemp/$rounded/unprocessed_beagle_input.tped";
		print $sys."\n";
		system($sys);
		$sys = "mv $dirtemp/$rounded/$array_arraytype[0]/unprocessed_beagle_input.tfam $dirtemp/$rounded/unprocessed_beagle_input.tfam"; 
		print $sys."\n";
		system($sys);
		
		$sys = "mv $dirtemp/$rounded/$qc[0] $dirtemp/$rounded/sample.qc";
		print $sys."\n";
		system($sys);
	}
	else
	{
		system("cat $dirtemp/$rounded/totalsnps |sort|uniq > $dirtemp/$rounded/totalsnps1");
		system("mv $dirtemp/$rounded/totalsnps1 $dirtemp/$rounded/totalsnps");
		system("cut -f1-4 -d ' ' $dirtemp/$rounded/totalsamples |sort|uniq > $dirtemp/$rounded/totalsamples1");
		system("mv $dirtemp/$rounded/totalsamples1 $dirtemp/$rounded/totalsamples");
		system("$dir/bin/plink/plink --tfile $dirtemp/$rounded/unprocessed_beagle_input_ori --extract $dirtemp/$rounded/totalsnps  --keep $dirtemp/$rounded/totalsamples --transpose --recode --out $dirtemp/$rounded/unprocessed_beagle_input_after_qc");
		system("$dir/bin/plink/plink --tfile $dirtemp/$rounded/unprocessed_beagle_input_after_qc --missing --out $dirtemp/$rounded/miss");
	 	for($i=0;$i<@array_arraytype;$i++)
        	{
			$buff = "BUFF$i";
			open $buff,"$dirtemp/$rounded/$array_arraytype[$i].qc" or die "no file $array_arraytype[$i].qc" ;
			$sample_qc_intial = <$buff>;
			while($line=<$buff>)
			{
				chomp($line);
				@array = ();
				@array = split(/\t/,$line);
				if(exists(${"Array_sampleqc"}{$array[0]}))
				{
					${"Array_sampleqc"}{$array[0]} = ${"Array_sampleqc"}{$array[0]}."\t".$line;
				}
				else
				{
					${"Array_sampleqc"}{$array[0]} = $line;		
				}
			}
			close($buff);
			open $buff,"$dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tfam " or die "no file $dirtemp/$rounded/$array_arraytype[$i]/unprocessed_beagle_input.tfam ";
			 while($line=<$buff>)
                        {
                                chomp($line);
                                @array = ();
                                @array = split(/ /,$line);
                                ${"Sex_sampleqc$i"}{$array[1]} = $array[4];
				#print "Sex_sampleqc$i\t".${"Sex_sampleqc$i"}{$array[1]}."\tSex_sampleqc$i"."\n";
                        }
                        close($buff);

		}
#		open(TFAM,"$dirtemp/$rounded/unprocessed_beagle_input_after_qc.tfam");
#		while($line = <TFAM>)
#		{
#			chomp($line);
#			@array = ();
#			@array  = split(" ",$line);
#			$tfam{$array[1]} =$line;
#		}
		close(TFAM);
		open(IMISS,"$dirtemp/$rounded/miss.imiss");
		$line = <IMISS>;
		while($line = <IMISS>)
		{	
			chomp($line);
			$line =~ s/\s+/ /g;
			@array = ();
			@array = split(" ",$line);
			$imiss{$array[1]} = $array[3]."\t".$array[4];
		}
		open(NEWSAMPLEQC,">$dirtemp/$rounded/sample.qc");
		print NEWSAMPLEQC $sample_qc_intial;
		while(($key,$value) =each %{"Array_sampleqc"})
		{
			chomp($key);
			chomp($value);
			@array = ();
			@array = split(/\t/,$value);
			  	
			if(@array_arraytype*8 != @array)
			{
				die "number of sample QC files not equal to number of platforms num_of_arraytypes".@array."\t".@array_arraytype."\n";
			}
			$newsampleqc = $array[0];
			$k =1;
			$sampleflag =$array[$k];
			$miss_pencent = $array[$k+1]/$array[$k+2];
			#$miss_pencent = $array[$k+2];
			$sex_ind = $array[$k+4]."\t".$array[$k+5]."\t".$array[$k+6];
			$sex_code{$key} = 0;			
			for($i=1;$i<@array_arraytype;$i++)
			{
				$k =$k+8;
				$sampleflag = $sampleflag."\t".$array[$k];
				$miss_pencent1 = $array[$k+1]/$array[$k+2];
				
				#$miss_pencent1 = $array[$k+2];
				$sex_ind1 = $array[$k+4]."\t".$array[$k+5]."\t".$array[$k+6];
				if($miss_pencent1 < $miss_pencent && $array[$k]  =~ m/Y/)
				{
					$miss_pencent = $miss_pencent1;
					$sex_ind = $sex_ind1;
					$sex_code{$key} = $i;
				}
			}	
			#print "sexcode".$sex_code{$key}."\n";
			if($sampleflag =~ m/Y/)
			{
				$newsampleqc =$newsampleqc."\tY\t".$imiss{$key}."\t"."1\t".$sex_ind."\n";
			}
			else
			{
				$newsampleqc =$newsampleqc."\tN\t".$imiss{$key}."\t"."0\t"."0\t"."NA\t"."NA\n";	
			}
			#print NEWSAMPLEQC $newsampleqc;	
			$hash_final_sampleqc{$key} = $newsampleqc;		
		}
		close(BUFF);
		open(BUFF,"$dirtemp/$rounded/$array_arraytype[0].qc");
		while($line = <BUFF>)
		{	chomp($line);
			@array =();
			@array = split("\t",$line);
			print NEWSAMPLEQC $hash_final_sampleqc{$array[0]};
		}
		open(NEWTFAM,">$dirtemp/$rounded/unprocessed_beagle_input_after_qc_final.tfam");
		open(BUFF,"$dirtemp/$rounded/unprocessed_beagle_input_after_qc.tfam");
		while($line = <BUFF>)
		{
			chomp($line);
			#print $line."\n";
			@array =();
 			@array = split(" ",$line);
			if(exists(${"Sex_sampleqc$sex_code{$array[1]}"}{$array[1]}))
			{
				$array[4] = ${"Sex_sampleqc$sex_code{$array[1]}"}{$array[1]};
			}
			else
			{

				die "not value exists Sex_sampleqc$sex_code{$array[1]} $array[1]\n";
			}
			chomp($array[4]);
			$array[4] =~ s/\s//g;
			$array[4] =~ s/\t//g;
			$array[4] =~ s/\n//g;
			$array[4] =~ s/\r//g;
			#print "Sex_sampleqc$sex_code{$array[1]}\t $array[1]"."\n";
			$line = join(" ",@array);
			print NEWTFAM $line."\n";
		}
		system("mv $dirtemp/$rounded/unprocessed_beagle_input_after_qc_final.tfam $dirtemp/$rounded/unprocessed_beagle_input.tfam");
		system("mv $dirtemp/$rounded/unprocessed_beagle_input_after_qc.tped $dirtemp/$rounded/unprocessed_beagle_input.tped");
	}	
#}



#forward strand mapping
@array_arraytype = split(',',$arraytype);
if(@array_arraytype > 1)
{
	$refpath = "$dir/bin/FwdStrand_Naresh/ref/".$array_arraytype[0]."/FWD_REF_INFO";
	for($ij=1;$ij<@array_arraytype;$ij++)
	{
		$refpath =$refpath.","."$dir/bin/FwdStrand_Naresh/ref/".$array_arraytype[$ij]."/FWD_REF_INFO";
	}
}
else
{
	$refpath = "$dir/bin/FwdStrand_Naresh/ref/$arraytype/FWD_REF_INFO";
}
#$cmd_system = "/bin/sh $dir/bin/FwdStrand_Naresh/bin/RunFwd.sh $dir/bin/FwdStrand_Naresh $dir/temp/$rounded/ $dir_or/$tped $dir/temp/$rounded/logfile $arraytype $refpath";
$cmd_system = "/bin/sh $dir/bin/FwdStrand_Naresh/bin/RunFwd.sh $dir/bin/FwdStrand_Naresh $dirtemp/$rounded/ $dirtemp/$rounded/unprocessed_beagle_input.tped $dirtemp/$rounded/logfile $arraytype $refpath $dir/bin/complete1KGSNPS_split";
print $cmd_system."\n";
system($cmd_system);
#die;
wait;
#system("cp $tfam  $dir/temp/$rounded");
#@tfam_ar=split(/\//,$tfam);
#$tfam = pop(@tfam_ar);
system("mv $dirtemp/$rounded/fwdStrandResults_1kg.tped $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped");
system("mv $dirtemp/$rounded/flippedInd_1kg.txt  $dirtemp/$rounded/fwdStrandFlippedInd_unprocessed_beagle_input.txt");
wait;
$FW_TP = "fwdStrandResults_$ori_tped";
$FW_TXT = "fwdStrandFlippedInd_$ori_tped";
$FW_TXT =~ s/tped/txt/g;
$SAMP = "sampleqc_$ori_tped";
$FW_TF = $FW_TP;
$FW_TF =~ s/tped/tfam/g;
$file = "$dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped";
if (!(-e $file)) {
        die "Fwd output tped file does not exist $file\n";
}
$file = "$dirtemp/$rounded/fwdStrandFlippedInd_unprocessed_beagle_input.txt";
if (!(-e $file)) {
        die "Fwd output ind file does not exist $file\n";
}
$file = "$dirtemp/$rounded/unprocessed_beagle_input.tfam";
if (!(-e $file)) {
        die "tfam output file does not exist $file\n";
}
$file = "$dirtemp/$rounded/sample.qc";
if (!(-e $file)) {
        die "sample qc does not exist $file\n";
}

#snps filetered out
system("perl $dir/bin/perl_generate_snplist_discarded_snpqc.pl $dir_or/$ori_tped $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped $dirtemp/$rounded/snps_failed_QC");
print "$dirtemp/$rounded/fwdStrandResults_$tped $dirtemp/$rounded/$FW_TXT $dirtemp/$rounded/unprocessed_beagle_input.tfam\n";
#system("mv $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped  $outdir/$FW_TP");
#system("mv $dirtemp/$rounded/fwdStrandFlippedInd_unprocessed_beagle_input.txt $outdir/$FW_TXT");
#system("mv $dirtemp/$rounded/unprocessed_beagle_input.tfam $outdir/$FW_TF");
#system("mv $dirtemp/$rounded/sample.qc $outdir/$SAMP");

$sys = "cp $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped  $outdir/$FW_TP";
print $sys."\n";
system($sys);
$sys = "mv $dirtemp/$rounded/fwdStrandFlippedInd_unprocessed_beagle_input.txt $outdir/$FW_TXT";
print $sys."\n";
system($sys);
$sys = "cp $dirtemp/$rounded/unprocessed_beagle_input.tfam $outdir/$FW_TF";
print $sys."\n";
system($sys);
$sys = "cp $dirtemp/$rounded/sample.qc $outdir/$SAMP";
print $sys."\n";
#system($sys);
system("mv $dirtemp/$rounded/snps_failed_QC $outdir/");
#system("rm -rf $dirtemp/$rounded/");


#include Structure

system("mv $dirtemp/$rounded/unprocessed_beagle_input.tfam $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tfam");

system("cut -f2 -d ' ' $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input.tped > $dirtemp/$rounded/snps_inputdata");
system("cut -f2 $dir/bin/reference.bim > $dirtemp/$rounded/distanct_snps");
system("cat $dirtemp/$rounded/snps_inputdata $dirtemp/$rounded/distanct_snps |sort|uniq -d > $dirtemp/$rounded/snps_common_input_and_distainct_selected"); 
system("$dir/bin/plink/plink --bfile $dir/bin/reference --keep $dir/bin/selected_1000genome_hapmap_refknown --extract $dirtemp/$rounded/snps_common_input_and_distainct_selected --make-bed --out $dirtemp/$rounded/reference");
system("$dir/bin/plink/plink --tfile $dirtemp/$rounded/fwdStrandResults_unprocessed_beagle_input --extract $dirtemp/$rounded/snps_common_input_and_distainct_selected --make-bed --out $dirtemp/$rounded/inputdata");

#fliping snp in input comparing to hapmap
open(BUFF,"$dirtemp/$rounded/inputdata.bim") or die "no file found $dirtemp/$rounded/inputdata.bim\n";
open(REF,"$dirtemp/$rounded/reference.bim") or die "no file found $dirtemp/$rounded/reference.bim\n";
open(WRBUFF1,">$dirtemp/$rounded/combined_ref_input.missnp") or die "not able to write $dirtemp/$rounded/combined_ref_input.missnp\n";
open(WRBUFF2,">$dirtemp/$rounded/combined_ref_input.drop") or die "not able to write $dirtemp/$rounded/combined_ref_input.drop\n";
$drop = 0;
$flip = 0;
$totalnum =0;
while(<REF>)
{
	chomp($_);
	@ref = split("\t",$_);
	$REF{$ref[1]} = $_;
	$totalnum++;
}
if($totalnum < 500)
{
	die "not enough snps present to do structure analysis\n";
}
while(<BUFF>)
{
	#$line = <REF>;
	chomp($_);
	#chomp($line);
	@input = split("\t",$_);
	if(exists($REF{$input[1]}))
	{
		@ref = split("\t",$REF{$input[1]});
	}
	else
	{
		die "snp in input but not in the hapmap reference : structure :dataset prep step\n";
	}
#if monomorphic condition	
	if($input[4] ne "0")
	{
		if(($input[5] eq $ref[5] && $input[4] eq $ref[4]) ||($input[4] eq $ref[5] && $input[5] eq $ref[4]))
		{
		}
		else
		{
			$input[5] =~ tr/ATGC/TACG/;
			$input[4] =~ tr/ATGC/TACG/;
			if(($input[5] eq $ref[5] && $input[4] eq $ref[4]) ||($input[4] eq $ref[5] && $input[5] eq $ref[4]))
			{
				print WRBUFF1 $input[1]."\n";
				$flip++;	
			}
			else
			{
				print WRBUFF2 $input[1]."\n";
				$drop++;
			}		
		}	
	}
	else
	{
		if($input[5] eq $ref[5] || $input[5] eq $ref[4])
		{
		}
		else
		{
			$input[5] =~ tr/ATGC/TACG/;
			if($input[5] eq $ref[5] || $input[5] eq $ref[4])
			{
				print WRBUFF1 $input[1]."\n";
				$flip++;	
			}
			else
			{
				print WRBUFF2 $input[1]."\n";
				$drop++;
			}		
		}	
	
	}
}	
close(WRBUFF1);
close(WRBUFF2);
close(REF);
close(BUFF);

if($flip >1 && $drop > 1)
{
	system("$dir/bin/plink/plink --bfile $dirtemp/$rounded/inputdata --exclude $dirtemp/$rounded/combined_ref_input.drop --flip $dirtemp/$rounded/combined_ref_input.missnp --make-bed --out $dirtemp/$rounded/inputdata1");
	system(" mv $dirtemp/$rounded/inputdata1.bim $dirtemp/$rounded/inputdata.bim");
	system(" mv $dirtemp/$rounded/inputdata1.bed $dirtemp/$rounded/inputdata.bed");
	system(" mv $dirtemp/$rounded/inputdata1.fam $dirtemp/$rounded/inputdata.fam");	
	system("$dir/bin/plink/plink --bfile  $dirtemp/$rounded/reference  --bmerge $dirtemp/$rounded/inputdata.bed $dirtemp/$rounded/inputdata.bim $dirtemp/$rounded/inputdata.fam --recode12 --output-missing-genotype -9 --out $dirtemp/$rounded/combined_ref_input"); 
}
elsif($flip >1)
{
	system("$dir/bin/plink/plink --bfile $dirtemp/$rounded/inputdata  --flip $dirtemp/$rounded/combined_ref_input.missnp --make-bed --out $dirtemp/$rounded/inputdata1");
	system(" mv $dirtemp/$rounded/inputdata1.bim $dirtemp/$rounded/inputdata.bim");
	system(" mv $dirtemp/$rounded/inputdata1.bed $dirtemp/$rounded/inputdata.bed");
	system(" mv $dirtemp/$rounded/inputdata1.fam $dirtemp/$rounded/inputdata.fam");	
	system("$dir/bin/plink/plink --bfile  $dirtemp/$rounded/reference  --bmerge $dirtemp/$rounded/inputdata.bed $dirtemp/$rounded/inputdata.bim $dirtemp/$rounded/inputdata.fam --recode12 --output-missing-genotype -9 --out $dirtemp/$rounded/combined_ref_input"); 

}
elsif($drop > 1)
{
	system("$dir/bin/plink/plink --bfile $dirtemp/$rounded/inputdata --exclude $dirtemp/$rounded/combined_ref_input.drop  --make-bed --out $dirtemp/$rounded/inputdata1");
	system(" mv $dirtemp/$rounded/inputdata1.bim $dirtemp/$rounded/inputdata.bim");
	system(" mv $dirtemp/$rounded/inputdata1.bed $dirtemp/$rounded/inputdata.bed");
	system(" mv $dirtemp/$rounded/inputdata1.fam $dirtemp/$rounded/inputdata.fam");	
	system("$dir/bin/plink/plink --bfile  $dirtemp/$rounded/reference  --bmerge $dirtemp/$rounded/inputdata.bed $dirtemp/$rounded/inputdata.bim $dirtemp/$rounded/inputdata.fam --recode12 --output-missing-genotype -9 --out $dirtemp/$rounded/combined_ref_input"); 

}
else
{
	system("$dir/bin/plink/plink --bfile  $dirtemp/$rounded/reference  --bmerge $dirtemp/$rounded/inputdata.bed $dirtemp/$rounded/inputdata.bim $dirtemp/$rounded/inputdata.fam --recode12 --output-missing-genotype -9 --out $dirtemp/$rounded/combined_ref_input"); 
}
#preparing the structure program
system("mkdir $dirtemp/$rounded/workdir");
system("cut -d ' ' -f2,6- $dirtemp/$rounded/combined_ref_input.ped >  $dirtemp/$rounded/workdir/temp.geno");

#renaming the population column(second) hapmap  samples in the structure dataset
open BUFF,"$dir/bin/referencesamples" or die "no file $dir/bin/referencesamples found\n";
while($line = <BUFF>)
{
        chomp($line);

        @array = split(/\t/,$line);
        #print "@array\n";
        $hash{$array[0]} = $array[2];
}

open BUFF,"$dirtemp/$rounded/workdir/temp.geno" or die "no file infile found\n";
open(WRBUFF,">$dirtemp/$rounded/workdir/temp1.geno") or die "unable to write file\n";
$total_samp_num = 0;
$totalnum = 0;
while($line = <BUFF>)
{
		if($total_samp_num == 0)
		{
			@totalnum = split(" ",$line);
			$totalnum = (@totalnum -2)/2;
		}
		$total_samp_num++;
        $line =~ s/\t+/ /g;
        $line =~ s/\s+/ /g;
        chomp($line);
        @array = split(/ /,$line);
        if(exists($hash{$array[0]}))
        {
                $array[1] = $hash{$array[0]};
        }
        else
        {
                $array[1] = -9;
        }
        $line = join(" ",@array);
        print WRBUFF $line."\n";
}
close(BUFF);
close(WRBUFF);
system("mv $dirtemp/$rounded/workdir/temp1.geno $dirtemp/$rounded/workdir/temp.geno");
#creating the extraparam file
system("cp $dir/bin/structure.extraparams $dirtemp/$rounded/workdir/");
#creating the main param file
open(WRBUFF,">$dirtemp/$rounded/workdir/structure.mainparams") or die "unable to write main param file\n";
print WRBUFF '#define INFILE'." $dirtemp/$rounded/workdir/temp.geno\n";
print WRBUFF '#define OUTFILE '."$dirtemp/$rounded/workdir/temp_structure.out\n";
print WRBUFF '#define NUMINDS '."$total_samp_num\n";
print WRBUFF '#define NUMLOCI '."$totalnum\n";
print WRBUFF '#define LABEL     1'."\n";
print WRBUFF '#define POPDATA   1'."\n";
print WRBUFF '#define POPFLAG   0'."\n";
print WRBUFF '#define PHENOTYPE 0'."\n";
print WRBUFF '#define EXTRACOLS 0'."\n";
print WRBUFF '#define PHASEINFO 0'."\n";
print WRBUFF '#define MISSING      -9'."\n";
print WRBUFF '#define PLOIDY       2'."\n";
print WRBUFF '#define ONEROWPERIND 1'."\n";
print WRBUFF '#define MARKERNAMES  0'."\n";
print WRBUFF '#define MAPDISTANCES 0'."\n";
print WRBUFF '#define MAXPOPS    2'."\n";
print WRBUFF '#define BURNIN    4000'."\n";
print WRBUFF '#define NUMREPS   4000'."\n";
system("$dir/bin/structure -m $dirtemp/$rounded/workdir/structure.mainparams -e $dirtemp/$rounded/workdir/structure.extraparams -K 3"); 

if(!(-e "$dirtemp/$rounded/workdir/temp_structure.out_f"))
{
	die "out file $dirtemp/$rounded/workdir/temp_structure.out_f does not exists\n";
}
$sys = "perl $dir/bin/perl_process_after_result_structure_out.pl $dirtemp/$rounded/workdir/temp_structure.out_f $dirtemp/$rounded/QC_output_structure $dir/bin/newreference_info";
print "$sys\n";
system($sys);
#combining result with QC result with struture result
open(BUFF,"$dirtemp/$rounded/QC_output_structure") or die " no file $dirtemp/$rounded/QC_output_structure structure output exists\n";
open(INPUT,"$dirtemp/$rounded/inputdata.fam") or die  " no file $dirtemp/$rounded/inputdata.fam exists\n";
open(QCIN,"$dirtemp/$rounded/sample.qc " ) or die  " no QC output file $dirtemp/$rounded/sample.qc exists\n";
open(QCOUT,">$dirtemp/$rounded/sample_out_strctureQC" ) or die  " no able to write structure output file\n";
$line=<BUFF>;
chomp($line);
$line =~ s/Sample\t//g;
$qcin = <QCIN>;
chomp($qcin);
print QCOUT "$qcin\t$line\n";
while(<BUFF>)
{
	chomp($_);
	@line = split(" ",$_);
	shift(@line);
	$_ = join(" ",@line);
	$input = <INPUT>;
	@input = split(" ",$input);
	$structure{$input[1]}= $_;
	#print "$input\t$_\n";
}

while($qcin = <QCIN>)
{
	chomp($qcin);
	@qcin = split("\t",$qcin);
	if(exists($structure{$qcin[0]}))
	{
		print QCOUT "$qcin\t".$structure{$qcin[0]}."\n";
	}
	else
	{
		print QCOUT "$qcin\tNA\tNA\tNA\n";
	}
}
close(BUFF);
close(INPUT);
close(QCIN);
close(QCOUT);

system("mv $dirtemp/$rounded/sample_out_strctureQC $outdir/$SAMP");
if($dirtemp =~ m/\w/ && $rounded =~ m/\w/ )
{
	#system("rm -rf $dirtemp/$rounded");
}
